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ABSTRACT 

The  effect  of  activation  energy  on  the  steady  state  solutions  exhibited 
by  the  nonadiabatic  tubular  reactor  is  investigated  by  applying  numerical 
bifurcation  techniques  to  the  model  equations.  As  the  activation  energy  is 
increased,  the  solution  branch  becomes  more  complex  and  finally  evolves  into  a 
multiplicity  pattern  with  regions  of  one,  three,  five  and  seven  solutions. 

Only  the  states  of  lowest  and  highest  temperature  are  found  to  be  stable. 

This  work  confirms  recent  results  using  large  activation  energy  asymptotics 
and  links  these  results  to  previous  numerical  studies  of  reactor  multiplicity 
employing  low  to  moderate  values  of  activation  energy. 
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SIGNIFICANCE  AND  EXPLANATION 


The  solutions  exhibited  by  chemical  reactor  models  not  only  predict 
performance  but  in  the  case  of  multiple  steady  states  also  elucidate  reactor 
dynamics  which  arise  due  to  transitions  between  stable  solutions.  We  trace 
the  response  curves  for  the  tubular  reactor  with  an  A  ♦  B  reaction  as  the 
activation  energy  y  (the  model  parameter  which  controls  the  effect  of 
temperature  on  reaction  rate)  is  increased  from  low  to  extremely  large 
values.  These  curves  become  quite  complex  as  y  is  assigned  large  values, 
however  the  dynamic  capabilities  of  reactor  do  not  increase  with  y  since 
only  the  solutions  of  highest  and  lowest  temperature  are  stable. 


The  responsibility  for  the  wording  and  views  expressed  in  this  descriptive 
summary  lies  with  MRC ,  and  not  with  the  author  of  this  report. 
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INTRODUCTION 

The  existence  of  multiple  steady  states  in  the  tubular  reactor  with  axial 
mixing  was  first  illustrated  by  the  numerical  computations  of  van  Heerden 

(I) .  Many  papers  have  followed  his  lead  by  theoretically  and  numerically 
investigating  the  uniqueness,  stability,  sensitivity,  and  structure  of  the 
solutions  exhibited  by  this  fundamental  type  of  reactor.  Theoretical 
examinations  have  centered  around  the  determination  of  criteria  for  uniqueness 
and  stability  in  terms  of  the  model  parameters.  Several  analytical  methods 
have  been  employed:  first-order  lumping,  Liaponov  functions,  comparison 
theorems,  and  fixed-point  analyses.  Their  application  has  shown  that  in 
general  the  solutions  exhibited  by  the  reactor  are  unique  for  sufficiently 
large  Peclet  numbers,  large  heat  transfer  coefficients,  small  reactor  length, 
and  small  Damkohler  number.  These  conclusions  are  given  in  the  reviews  of 
Varma  and  Aris  (2)  and  Schmitz  (3)  and  in  the  works  of  Hlavacek  and  Hofmann 
(4),  Luss  et  al.  (5-7),  Varma  and  Amundson  (8-10)  and  McGowin  and  Perlmutter 

(II) . 

Multiple  solutions  have  been  numerically  studied  by  applying  variations 
of  the  shooting  and  collocation  methods  to  the  reactor  equations  using  model 
parameters  well  outside  the  range  given  by  the  uniqueness  criteria  in  the 
above  references.  The  investigations  of  Varma  and  Amundson  (12-13),  Hlavacek 
(14-16),  McGowin  and  Perlmutter  (17),  and  Jensen  and  Ray  (18)  have  determined 


♦Department  of  Mathematics,  Colorado  State  University 
Fort  Collins,  Colorado  80623 


Sponsored  by  the  United  States  Army  under  Contract  No.  DAAG2O-80-C-004 1 . 


the  existence  and  stability  of  one,  three,  and  five  steady  state  solutions  and 
have  examined  the  effect  of  the  system  parameters  on  the  range  of 
multiplicity.  Although  sustained  oscillations  were  also  demonstrated  in  the 
nonadiabatic  case,  the  above  studies  were  primarily  directed  toward  the 
multiplicity  problem  and  the  associated  hysteresis  and  jump  phenomena  which 
occur  at  transitions  between  stable  operating  states  in  the  reactor.  Tapila 
and  Poore  (19)  have  recently  classified  the  steady  state  solutions  for  all 
system  parameters  in  the  asymptotic  limit  of  large  activation  energy.  These 
perturbational  methods  show  that  the  reactor  may  in  fact  exhibit  seven  steady 
states  for  higher  activation  energies. 

The  goal  of  this  work  is  to  show  how  the  one,  three  and  five  steady 
states,  which  were  studied  extensively  in  the  seventies,  evolve  into  the  newly 
discovered  seven  steady  states.  The  present  study  then  confirms  the  work  of 
Kapila  and  Poore  and  ties  it  to  the  work  of  the  previous  investigators. 

MATHEMATICAL  MODEL  AND  NUMERICAL  TECHNIQUES 

The  steady  state  conservation  equations  describing  an  A  B  reaction  in 
the  nonadiabatic  tubular  reactor  with  axial  mixing  are  written  below  in 
dimensionless  form 
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These  equations  are  subject  to  the  following  boundary  conditions 
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The  variables  c  and  T  represent  the  dimensionless  concentration  of 

reactant  A  and  reactor  temperature.  The  parameters  Pe  ,Pe.  ,&,T  ,B,D,  and 

m  h  c 

Y  denote  respectively  the  Peclet  number  for  mass  transfer,  the  Peclet  number 
for  heat  transfer,  the  heat  transfer  coefficient,  the  coolant  temperature,  the 
heat  of  reaction,  the  Damkohler  number  and  the  activation  energy. 

We  investigate  reactor  multiplicity  by  computing  the  steady  state 
solutions  of  the  above  model  as  the  Damkohler  number  is  varied.  The  boundary 
value  problem  is  discretized  by  a  fourth-order  finite  difference  scheme  due  to 
Stepleman  (20),  and  the  resulting  nonlinear  algebraic  equations  are  solved  by 
Keller's  arclength  modification  of  the  Euler-Newton  continuation  procedure 
(21).  Since  these  methods  are  discussed  in  the  literature  (22,  23),  we  forego 
any  further  discussion  of  these  techniques.  The  stability  of  the  steady 
states  is  determined  by  computing  the  eigenvalues  of  the  discretized  system 
via  the  QZ  algorithm. 


MTTLTIPLICITY  PATTERNS 

The  results  of  our  computations  are  represented  on  response  curves  which 

represent  the  maximum  reactor  temperature  (the  largest  value  of  T  found  on 

the  computed  temperature  profile)  plotted  against  the  Damkohler  number.  We 

examine  the  effect  of  the  activation  energy  on  the  reactor  multiplicity  by 

tracing  the  response  curve  as  y  varies  from  25  to  225  while  the  remaining 

parameters  are  fixed  at  Pe  =  Pe  =  1.0,  6=  4.0,  B  =  0.5  and  T  = 

m  n  c 


1.0 


This  set  of  parameters  falls  into  the  range  for  which  seven  steady  states  are 

predicted  by  the  analysis  of  Kapila  and  Poore.  More  specifically,  let 

r  ■  ^  Pe.  [1  ±  (1  +  46/PeJ/2]  and 
I  2  n  n 

a( 0, Pe,)  =  (exp(r.)  -  exp(r  ))/(r  exp(r)  -  r  exp(r  )).  If  6  denotes  the 
n  +  -  +  +  -  -  c 

value  of  B  for  which  Ba(  B/Pe,  )  =  1,  then  in  the  asymptotic  limit  of  large 

h 

activation  energy,  there  exists  for  some  range  of  the  Damkohler  number,  three 

solutions  for  $  =  0,  five  solutions  for  0  <  B  <  8  and  seven  solutions  for 

c 

6  >  6  (19).  (For  Pe.  =  0,  0  is  approximately  1.4  while  for  Pe.  larger 

c  h  c  n 

than  4,  0^  is  slightly  larger  than  1  +  Pe^.)  Although  we  present  the 

transition  to  seven  steady  states  for  only  one  set  of  parameters,  the  results 
are  quite  similar  for  other  sets  satisfying  the  above  criteria. 

In  Figure  1,  we  show  that  the  reactor  model  yields  a  branch  of  unique, 
stable  steady  state  solutions  for  y  =  25.  As  y  is  increased  to  40,  a 
region  of  three  solutions  appears  for  a  very  narrow  range  of  the  Damkohler 
number.  A  close  examination  of  this  solution  branch  indicates  that  all  three 
solutions  may  be  unstable  since  exchanges  in  stability  do  not  occur  at  the 
static  bifurcation  points  but  rather  at  the  Hopf  bifurcation  points  where 
periodic  solutions  bifurcate  from  the  solution  branch  and  surround  the 
unstable  steady  states  (24,  19). 

The  remaining  solution  branches  in  Figure  1  all  exhibit  five  solutions 
over  some  range  of  D.  An  inner  loop,  which  grows  as  y  is  increased, 
appears  on  the  intermediate  section  of  the  usual  S-shaped  curve,  and  this  loop 
gives  rise  to  the  five  steady  states.  For  each  of  these  solution  branches  the 
first  exchange  of  stability  occurs  at  the  lower  limit  point,  and  all  of  the 
intermediate  solutions  are  unstable  until  a  Hopf  bifurcation  point  is 
encountered  above  and  to  the  right  of  the  upper  limit  point.  We  also  note 
that  the  y  =  75  branch  has  the  identical  shape  of  the  multiplicity  pattern 
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predicted  by  Kapila  and  Poore  (19)  for  0  <  8  <  8C«  (Our  computations  with 

Pe  =  1.0,  g  =  1.0,  B  =  0.5  and  y  =  100  confirm  the  Kapila  and  Poore 

predictions  for  0  <  g  <  8  but  are  not  presented  here  because  of  their 

c 

similarity  with  y  =  75  results  of  Figure  1.) 

In  Figure  2,  y  is  increased  to  100  and  125,  and  seven  solutions  appear 
in  both  cases  when  two  limit  points  arise  along  the  lower,  right-hand  portion 
of  the  inner  loop.  This  region  is  barely  detectable  for  y  *  100  and  is 
slightly  broadened  for  y  =  125  as  this  S  shaped  perturbation  of  the  inner 
loop  is  more  clearly  defined.  The  lower  two  solutions  on  this  perturbation 
correspond  to  the  additional  solutions  established  by  Kapila  and  Poore.  The 
stability  results  for  the  branches  in  Figure  2  are  the  same  as  those  for 
y  =  75,  where  all  the  solutions  except  the  lowest  and  highest  temperature 
states  are  unstable. 

The  exclusion  of  the  lower  limit  points  (D  =  0.042  for  y  =  100  and 

D  =  0.032  for  y  =  125)  from  Figure  2  indicates  the  effect  of  large  y  on 

the  location  of  the  multiplicity  regions.  It  is  clear  from  the  definition 

of  D  (D  is  proportional  to  e  (3).)  that  increasing  y  to  such  large 

values  must  shift  the  solution  branch  rather  dramatically  to  the  left. 

The  result  of  increasing  y  to  225  is  illustrated  in  Figure  3  where 

seven  solutions  are  shown  for  a  narrow  range  of  D 
-10  -9 

(7.0  x  10  <  D  <  1.35  x  10  )  and  five  solutions  are  shown  over  a 

-12  -6 

relatively  wide  range  (5.0  x  10  <  D  <  5.5  x  10  )  as  the  perturbation  on 

the  inner  loop  becomes  quite  exaggerated.  The  maximum  temperature  of  the 
newly  discovered  solutions  becomes  identical  along  the  lower  portion  of  the 
perturbation,  but  the  locations  of  the  maximums  occur  at  different  points 
along  the  reactor.  The  solution  branch  for  this  extremely  large  value  of  y 
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DAMKOHLER  NUMBER 


shows  excellent  agreement  with  the  asymptotic  results  for  large  activation 
energy  presented  by  Kapila  and  Poore  [19]. 

CONCLUSIONS 

We  have  illustrated  the  effects  of  activation  energy  on  the  steady  stat 
multiplicity  of  the  tubular  reactor.  As  the  activation  energy  increases,  th 
solution  branch  becomes  more  complex  am1  evolves  into  a  pattern  of  seven 
solutions  of  which  only  two  (states  of  lowest  and  highest  temperature)  are 
stable . 
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